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ABSTRACT 

We present new interferometric data obtained with MIDI (MID infrared Interferometric instru- 
ment) for the Seyfert II galaxy NGC 1068, with an extensive coverage of sixteen uv points. 
These observations resolve the nuclear mid-infrared emission from NGC 1068 in unprece- 
dented detail with a maximum resolution of 7 mas. For the first time, sufficient uv points have 
been obtained, allowing us to generate an image of the source using maximum entropy image 
reconstruction. The features of the image are similar to those obtained by modelling. We find 
that the mid-infrared emission can be represented by two components, each with a Gaussian 
brightness distribution. The first, identified as the inner funnel of the obscuring torus, is hot 
(-800K), 1.35 parsec long, and 0.45 parsec thick in FWHM at a PA=-42° (from north to 
east). It has an absorption profile different than standard interstellar dust and with evidence 
for dumpiness. The second component is 3 x 4 pc in FWHM with T=~300K, and we identify 
it with the cooler body of the torus. The compact component is tilted by ~ 45° with respect to 
the radio jet and has similar size and orientation to the observed water maser distribution. We 
show how the dust distribution relates to other observables within a few parsecs of the core of 
the galaxy such as the nuclear masers, the radio jet, and the ionization cone. We compare our 
findings to a similar study of the Circinus galaxy and other relevant studies. Our findings shed 
new light on the relation between the different parsec-scale components in NGC 1068 and the 
obscuring torus. 

Key words: techniques: interferometric -galaxies: Seyfert -galaxies(individual):NGC 1068- 
infrared: galaxies 



1 INTRODUCTION 

The A GN unification model ( AntonuccJ 1 19931 : lUrrv & Padovanil 
1 19951) explains the difference between type II Seyferts, which show 
only narrow emission lines, and type I, which show both narrow 
and broad emission lines, by stipulating that a torus-like structure 
surrounds the central engine and accretion disk. Thus, when ori- 
ented edge-on, the torus blocks the subparsec-sized broad emis- 
sion line region making the object appear as a type II, while har- 
bouring an unseen type I nucleus. Symmetry and angular momen- 
tum considerations suggest that the torus and the accretion disk are 
perpendicular to the bipolar jet, whose orientation is likely to be 
coupled closely to that of the accreted matter. The relative num- 
ber of Seyfert I/II dem ands that the torus be geometrically thick 
dOsterb rock & Shaw 1988), although thick structures orbiting com- 
pact objects will quickly lose their height and collapse into a thin 
disk. In order to maintain its stable torus structure, the random ve- 
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locities of the dust clouds need to be of an order similar to the or- 
bital speed, or a few hundred km/s. Since colliding dust clouds are 
destroyed at relative speeds as low as a few m eters per second, the 
dust clouds are not expected to last very long jKrolik & Beaelman 
1988). To overcome the difficulty of maintaining its inflated state 
different dust configuration ot her than a torus have b een introduced, 
such as the warped disk of ISanders et al.l il989l) . Alternatively, 
the dust has been descr ibed as a hydrodynamically driven outflow 
jKartie&Koniglll 19961) . 

NGC 1068 is considered to be the prototype Seyfert II galaxy, 
where the central source is obscured by dust. Its relatively small 
distance of 14.4 Mpc and high mid-infrared flux make the object 
ideally suited for the study of the nucleus and obscuring dust. Pre- 
vious MIDI observations of NGC 1068 revealed warm (320 K) sil- 
icate dust in a structure 2.1 parsecs thick and 3.4 parsecs in di- 
ameter, surrounding a small, hot (800K) compo nent whose shape 
and orientation could not be determined in detail Jjaffe et alj2004h 
(hereafter J04). These observations were, together with those of 
IWittkowski et al.l d20041) (W04), the first to spatially resolve direct 
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emission from the putative torus, and the first to show that a torus- 
like structure is indeed present in NGC 1068. 

The mid-infrared emission from the central region of 
NGC 1068, as seen by the largest single dish telescopes, is com- 
posed of an unresolved core plus extended emission. Deconvolved 
maps at 12/jm taken with VISIR (Imager a nd Spectromete r in the 
InfraRed at the Very Large Telescope) by iGalliano et alj d2005h 
reveal a set of discrete mid-infrared sources: seven in the north- 
eastern quadrant and five in the south-eastern qu adrant. The centra l 
source, observed with the 10m Keck telescope jBock et al.ll2OO0h . 
is extended by ~1" in the north-south direction and is unresolved 
in the east-west direction. About 2/3 of its flux can be ascribed to 
a core structure which is itself el ongated north-sou th in a tongue- 
shaped structure and according to lBock et al.1 ( feOOOh does not show 
a distinct unresolved core. Rec ent 12.8^m speckle images taken 
with VISIR in BURST mode bv lPoncelet et all J2007I) identify two 
major sources of emission at 12.8yum: a compact source (< 85 mas) 
and an elliptical source of size (< 140) mas x 1187 mas and 
PA~ 4°. It is the unresolved compact source which is associated 
with the dusty torus and is the subject of this study. Other parsec 
scale components in NGC 1068 that are related to the torus are 
the compact H 2 masers , appearing on the sky as a set of lin- 
ear spots (PA=-45°) that are misaligned with the jet and span a 
velocity range of 600 km/s with a sub-Keplerian velocity profile 
(v oc r~ 03 ). The kinematics of the maser spots indicate that the 
masers are located in a rotating disk with inner radiu s of -0.65 
pc and outer radius of ~1.1 pc (Gallimore et al. 2001). This disk 
traces the outer, colder part of the accretion disk where conditions 
allow for the formation of water masers. It was proposed earlier t hat 
these masers trace warm molecular dust ( Claus sen & Lo|[l~986l) . a 
conclusion supported by the observations presented here which re- 
veal the dusty torus to be at the outer edge of the maser disk. VLBA 
5 and 8.4GHz radio continuum images show a parsec-sized struc- 
ture with a major axis at PA=:-75°, most likely ind icating free-free 
emiss ion from hot (T=10 4 - 10 5 K) ionized gas ( Gallimo re et al.l 
120041) . According to the unification model, toroidal obscuration by 
dust is also responsible for the conical shape of the narrow-line re- 
gion. The ionization cone in NGC 1068 as seen in HST images is 
centred around PA=10°, while modelling of HST spectra, based on 
the kinematics of the gas indicate the ioni zation cone with an open- 
ing angle of 80° centred around PA=30° jPas et al.f2 006). roughly 
perpendicular to the maser spots. 



2 OBSERVATIONS, UV COVERAGE AND DATA 
REDUCTION 

MIDI is the mid-infrared interferometer located on Cerro Paranal in 
northern Chile and operated by the European Southern Observatory 
(ESO). MIDI functions as a classical Michelson interferometer, 
combining the light from two 8-meter un it telescopes (UTs). For a 
detailed description of the instrument see lLeinert et al .1 J2003I) . The 
main data product of MIDI is the correlated flux, which can be ex- 
plained as the spectrum of the source at a certain spatial resolution 
projected perpendicular to the baseline, that is to the projected sep- 
aration of the two telescopes. A total of 16 interferometric obser- 
vations were taken together with a single-dish total flux spectrum 
for each one. Table [T] contains the observation log and instrument 
settings. A map of the uv coverage is shown in Fig. [T] Due to the 
near zero declination of NGC 1068 the tracks are parallel to the u- 
axis, causing each observation to differ both by baseline length and 
position angle from its neighbours. The distribution of the tracks is 
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Figure 1. uv coverage [m] for NGC1068 and colour coded by date. Due to 
the ~0.0 declination of NGC1068, the uv tracks are parallel to the u-axis. 
uv coords [u,v] are complex conjugates of uv coords [-u,-v], and both are 
plotted since they are indistinguishable. 



Table 1. Log of the observations dates, unit telescopes used and baseline 
configuration in polar (D, PA) coordinates projected on sky. All observa- 
tions were done in GRISM mode with a spectral resolution of A/AA = 230. 
The same calibrator, HD10380, was used for all observations with identi- 
cal settings. Chopping with / = 1Hz, and throw amplitude of 15" at an 
angle a = 0° was applied during photometric measurement and acquisi- 
tion. Fringes were tracked in "olfset tracking" mode. The spatial resolution 
specified is calculated by 8 = A/(2.2B) where B is the projected baseline 
separation and A = Wfun. 





Date 


UT's 


D 


PA 


resolution 








[m] 


deg 


[mas] 


#1 


Nov. 11,2005 


UT1-UT4 


90.8 


50.0 


10.3 


#2 






109.2 


57.6 


8.6 


#3 






122. 


61.4 


7.7 


#4 






129.6 


63.2 


7.2 


#5 


Nov. 12, 2005 


UT1-UT3 


79.9 


10.0 


11.7 


#6 






84.7 


21.7 


11 


#7 






90.8 


29.9 


10.3 


#8 






96.8 


35.6 


9.7 


#9 


Nov. 13,2005 


UT2-UT3 


33.6 


15.5 


27.9 


#10 






36.6 


27.5 


25.6 


#11 






40.15 


36.1 


23.4 


#12 


Nov. 14, 2005 


UT2-UT4 


69.0 


79.8 


13.6 


#13 






89.98 


82.1 


10.4 


#14 






87.76 


82.0 


10.7 


#15 


July 12, 2006 


UT3-UT4 


52.3 


-67.3 


17.9 


#16 


Aug. 10, 2006 


UT2-UT4 


59.82 


78.2 


15.6 



concentrated in the second/fourth quarter of the uv plane, with only 
one observation on a perpendicular baseline, reflecting the fixed 
physical placement of the VLTI unit telescopes. After each inter- 
ferometric observation, a total flux measurement (i.e. single-dish 
spectrum) is taken from each telescope independently, and used to 
determine the visibility, i.e. the correlated interferometric flux di- 
vided by the geometric mean of the single dish fluxes. Only one 
such observation was rejected due to bad seeing, and the remaining 
15 fluxes were averaged together (Fig.[2]>. A typical MIDI observa- 



Resolving the obscuring torus in NGC 1068 with the power of infrared interferometry 3 



tion is recorded in 8000 frames, with an integration time of 0.036 
seconds for each frame. We examined the total flux recorded on 
the detector per frame, and rejected frames where the flux was un- 
stable. Data reduction was done by coherent visibility estimation 
wit h the Expert Work Station (EWS) package written in Leiden. 
See |jaffel (2004) for a detailed description of the coherent visibility 
estimation method. 



2.1 Calibration 

Calibration of the interferometric and photometric data was done 
using HD 10380 as the calibrator star for all observations. Each 
night's observation started with HD10380, and it was observed be- 
tween observations of NGC 1068. Thus, we obtained a calibration 
measurement for each uv point of NGC 1068. The calibrator's data 
reduction was done using the same parameters as for NGC 1068 
and went through the same tests for consistency as for the sci- 
ence target. Only one calibrator observation was rejected due to 
bad seeing. Calibration measurements taken just before and after 
NGC 1068 were compared and, if different, were averaged. For an 
unresolved point source, the correlated flux should equal the total 
flux. Assuming this, we can calibrate the correlated fluxes using the 
calibrators, without a need to use the photometry measurements, 
which are less reliable due to the strong atmospheric background. 
The calibrated fluxes are computed by dividing the correlated fluxes 
of the target by those of the calibrator, and then multiplying by the 
known flux of the calibrator. F or HD 10380 we used the spectral 
template of ICohen et al.l J 1999h . The calibrated photometric fluxes 
were produced by the same procedure for the photometric data. 



2.2 Correlated flux. vs. visibility 

Optical/infrared interferometric data is often presented as visibil- 
ity data, a dimensionless quantity, which measures the fraction of 
flux that is resolved for a given baseline. Here we prefer to use the 
correlated fluxes, the actual flux that the interferometer measures 
for a given baseline, as is common practice in radio interferometry. 
Determining the visibility, V, directly from the contrast of the in- 
terferometric fringes ( V = (/,„ a , - l m i„)l(l ma x + Imm)) is not practical 
since the amount of flux reaching the beam combiner from the two 
telescopes is generally not equal and therefore affects the contrast 
of the fringes. Given an isolated source, the visibilities equal the 
correlated fluxes divided by the total flux of the source. If, however, 
the source is embedded in emission on larger scales (as seems to be 
the case for NGC 1068), then the total flux determined by a larger 
field of view is contaminated by the large scale emission making 
it difficult to interpret. The correlated fluxes also free us from han- 
dling the strong atmospheric background, which is not correlated 
and is removed easily in the data reduction process. 



3 RESULTS 
3.1 The total flux 

The total flux of NGC 1068 is shown in Fig. [2] as measured by 
the average of two 8-meter UTs, each with a beam size of 250 mas 
. It is in rough agreement with previously published spectra from 
8m class telescopes. The only spectral feature seen in the total flux 
is the 9.7yum silicate feature, in absorption . We do not detect the 
very weak [SIV] emission line reported bv lRhee & LarkirJ d2006l) 
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Figure 2. Total (single-dish) flux. All 15 individual spectra are in grey with 
the mean given in black. The dash-dotted line is the model fit described in 
section |4~T1 which is composed of two flux components: component 1 (red- 
dashed) and component 2 (blue-dotted). Note the atmospheric O3 feature at 
9.7//m. 




Figure 3. Correlated fluxes #1-8 (solid black) and the model fit from section 
14.11 (red-dashed). The panels are sorted according to baseline length and 
the unit telescopes used, and are numbered according to the number in the 
observation log (TablefT}. 



an d Mason et al. (2006) with similar apertures. Space based instru- 
ment s show PAH features at A = 7.6, 9 and 10.9yum 1 Stu rm et al.l 
2000). These features are found to occur at larger distances from 
the nucleus. 
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Figure 4. Correlated fluxes #9-16 (solid black) and the model fit described 
in section PTTl (red-dashed). The blue filled curves show the contribution 
of component 2 to the correlated fluxes. The 'bump' at 9.7/im visible on 
fluxes #9-14,16 is due to the atmospheric O3 feature. The panels are sorted 
according to baseline length and the unit telescopes used, and are numbered 
according to the number in the observation log (Table [T}. 

3.2 The correlated (interferometric) fluxes 

The correlated fluxes, shown in Figures [3]S{4] display large varia- 
tions in the shape and depth of the silicate feature centred at 9.7//m, 
the only spectral feature present in the data. No traces of PAH 
molecules or the 12.9yum Ne line are seen. The clear relation seen in 
the data between the baseline length and the correlated flux shows 
that the source is resolved on every baseline. The 'bump' around 
9.7/im seen in baselines #9-14,16 (Figures [3]£(4]) is due to the at- 
mospheric ozone feature, centred 9.7//m. 

3.3 Discussion of the results 

The correlated fluxes in fact sample the SED of NGC 1068 about a 
range of different orientations and spatial resolutions. The struc- 
ture of the spectra, and in particular the 9.1 [im silicate feature, 
varies considerably with baselines, and is quite different from the 
total flux. The striking differences between the correlated and to- 
tal fluxes demonstrate how inadequate the total flux is as the sole 
means to constrain torus models. Indeed, numerous authors have 
attempted to fit tori models to the SED of NGC 1068, and many 
di fferent models fit the SED equally well, as was first demonstrated 
bv lGalliano et al.U2003l) . 

The depth of the silicate feature in the correlated fluxes is 
about 0.6 of the continuum level (as defined by the average of the 
flux at 8 and 13yum) and is independent of baseline length. This 
suggests that the emitting medium is compact and is located in a 



region much smaller than the absorbing medium. For the total flux, 
the silicate feature depth is smaller, ~ 0.3. The total flux can be 
successfully fitted with a combination of two grey body spectra, a 
warm (~300K) component, plus a hot (~800K) compact compo- 
nent, both behind uniform absorption screens as discussed in 34.11 
In general we find the total flux to be easily reproducible by our 
models described in 34. ll for a wide range of parameters. 

A note on phases and symmetry 

Apart from the correlated and total fluxes, one can also recover 
the differential phase (fdiffik) of the source from MIDI data. The 
differential phase is related to the source's true phase, cf>{k), by 
4>diff(k) = <p{k) - C x k, where k = 2n/A is the wave number and C 
is an unknown constant. Thus, any linear dependence of on k is 
removed from the differential phase in the data reduction process, 
and cannot be recovered. This missing information in the differen- 
tial phase makes it very challenging to include it in our models. 
We have therefore postponed treatment of the differential phase in 
NGC 1068 to a second paper. Without the phase information, only 
models or images which possess inversion symmetry can be ap- 
plied to the data. Throughout this paper, we assume a symmetric 
flux distribution. This is an instrumental limitation and we do not 
claim that the true source brightness possesses such symmetry. 



4 MODELLING 

As interferometric data is obtained in Fourier space, interpreting the 
data is not as straightforward as in direct imaging or spectroscopy. 
In infrared interferometry, where only a small number of the abso- 
lute value of the Fourier components can be measured, the results 
are often model dependent and difficult to interpret unambiguously. 
The effective resolution for each uv point is determined by A/B, 
the projected separation between the telescopes expressed in num- 
ber of wavelengths. This relation causes the spectral features of the 
source, which strongly depend on A, to become entangled with the 
spatial flux distribution as observed by telescopes of separation B. 
We have attempted to disentangle the two effects by using two dif- 
ferent modelling approaches. First, we treat the flux distribution as 
coming from two Gaussian components, each with a fixed size and 
orientation. Each component is assumed to be a grey body with a 
fixed temperature and situated behind a uniform absorption screen 
( 34. U . The model results show us that at 8/um only one component 
contributes to the correlated flux. Therefore, we chose to attempt to 
reconstruct an image at this wavelength, which is shown and dis- 
cussed in 34.21 The image agrees well with the modelling results. 
In our third approach ( 34.3t we look at each wavelength indepen- 
dently and fit a single Gaussian component to each, assuming no 
relation between the different wavelengths in the data, and assum- 
ing all the observed structures in the correlated fluxes arise from the 
different source sizes and orientation at each wavelength (channel 
by channel fitting). A good understanding of the observed proper- 
ties of NGC 1068 can then be achieved by combining the results 
from the different procedures. 

4.1 Two grey body model 

The first approach to disentangling structure and spectrum is to 
assume that changes to the correlated fluxes with wavelength are 
solely due to spectral effects, i.e. a grey body emission undergo- 
ing absorption while keeping the spatial parameters constant with 
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Figure 5. Compari son between the two components in our model and the 
12.5//m image of iBock et ajj <200Ch . taken with the 10m Keck telescope. 
The two components are plotted symmetrically, yet their relative position 
cannot be determined from the correlated fluxes. 



wavelength. In this model we treat the infrared emission as coming 
from two Gaussian grey body components of a fixed size and ori- 
entation, each one behind a uniform absorption screen. Two is the 
minimal number of components necessary to account well for the 
data, while the addition of more components does not improve our 
fit. The correlated flux of each component is given by: 



F C0TT {A,u,v) = j]BB(A)V(u/A,v/A)e~ 



(1) 



where BBG) is the emission from a black body of temperature T; V 
is the visibility of a Gaussian component whose major axis, minor 
axis, and PA are to be fitted; abs is the absorption curve of a chosen 
mineral or a combination of a few minerals as described below. 
The coefficient 77 is the grey body scaling factor and has a value 
0< rj <1, independent of wavelength. After computing the Fourier 
amplitudes corresponding to each uv point in the model, the two 
sets of amplitudes (one for each component) are then combined to 
produce the final correlated flux to be compared to the data: 



F corr \Fcorr\ £ F corr ^\ 



(2) 



The parameters /, m are for the two components to be moved with 
respect to each other in the plane of the sky, parallel to the RA 
and DEC directions, respectively. In practice we find that the data 
does not provide constraints for the relative positions /, m (see also 
94.1.11 1 and we have set them to zero (concentric components). 
Each component in the model has six parameters, making in twelve 
in total: the major and minor FWHM, the position angle 0, the BB 
temperature T, the optical depth t, and the scaling factor rj. The 
main features of the data can be reconstructed from values of the 
correlated flux at about six wavelengths, in the sense that if the 
model fits the data well at those wavelengths it will fit the rest well. 
We can then estimate that these 12 parameters are in fact fitting 
6 x 16 = 96 independent data points. The total flux is treated here 
as an extra uv point to be fitted, and is computed by setting the 
u and v coordinates in equation [2] to zero for each component. In 
practice we find that the total flux is easily reproduced in many of 
our models and for a different range of parameters, and may serve 
as an upper limit to the flux and size of each component. The most 
difficult part in our model is to determine the absorption curves in 
equationQ] We have selected several dust absorption templates, in- 
cluding stan dard galactic dust as observed towards the centre of the 
Milky Way temper et all |2004|) . as well as Ca 2 Al 2 Si0 7 , which 
was fo und to best fit the previous MIDI data (J04, Speck et al. 
(2000)). Our aim here is not to unambiguously determine the chem- 



o 
U 

g 0.6 h 



< 



1 1 1 1 1 1 1 1 1 


Al.,0, 




_ Ca>l 2 Si0 7 




— standard galactic 




-- fitted composite 


/ '/ 




1,1,1,1, 





8 9 10 , , , 11 12 13 

Alum] 

Figure 6. A comparison between the absorption template used in fitting 
the data (dashed) and the other spectral templates for dust particles com- 
monly found in astrophysical environments. The latter includes that of dust 
seen toward the centre of the Milky Way. These templates only include line 
absorption and do not include the continuum absorption. 



ical composition of the dust, which is not possible to do from the 
data, but rather to characterise the spectral template which will fit 
the data best, and to compare it with other known dust mineral 
templates. In this respect the template is any simple function f(X) 
which, when inserted into equation Q] will provide the best fit for 
the data. To achieve this, we first fit a combination of the minerals 
listed above and shown in Figure [6] The best-fit combination was 
improved on by allowing the value of r at several key wavelengths 
(i.e. 8, 9, 9.7, 1 1.5 and 13yum) to vary and then interpolated between 
the new values using cubic spline. This method is able to mimic a 
typical absorption template. The resulting template, which we des- 
ignate 'fitted composite', along with several other dust templates 
are plotted in Fig. [6] while the model parameters are given in Ta- 
ble. [2] We label the two components in our model as components 1 
and 2. Figures [3] and |4]plot the best fit model against the correlated 
fluxes, while Figure ff] plots the model's fit to the total flux plus 
the contribution from the individual components. The fitting was 
done by finding the least-squares solution using the Levenberg- 
Marquardt technique. The general trend is for the quality of the 
fit to become worse with increasing baseline length, as one might 
expect from the lower spatial resolution of a short-baseline obser- 
vation. In addition, component 2 only contributes to the correlated 
fluxes taken with shorter baselines, effectively increasing the num- 
ber of fit parameters. Some of the features in the correlated fluxes, 
in particular the 'step' seen most clearly in correlated fluxes #4 
and 5 at 10.5yum, cannot be reproduced by the simple model we 
propose. Our model will always generate smooth fluxes due to the 
smooth absorption templates we use, which not contain any such 
'steps'. Furthermore, introducing more complicated templates will 
not help to fit these unsmooth features since any change made to 
the templates will affect each correlated flux, while these unsmooth 
features only appear in some of the correlated fluxes. 

Seven out of sixteen modelled correlated fluxes 
(#7,8,9,10,11,15,16) are in excellent agreement with the data 
while the rest of the modelled fluxes are mostly in good to perfect 
agreement with the data either for A <10/itn or for A >10yum. The 
two components are illustrated in Figure [5] Given the range of 
spatial resolutions of the different baselines (spanning a factor ~4 
in length), the complex behaviour of the correlated fluxes and the 
striking simplicity of the model, it provides a surprisingly good fit 
to the data. 
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Table 2. Parameter fits for the two grey-body model of j|4. 1| The values 
here are based on the results of different models that differ in the chemical 
composition of the dust. To determine the errors we fixed the value of all but 
one parameter at a time and noted the range of values resulting in significant 
changes to the model. 

"The temperature and scaling factor, are not independent. Temperatures 
as high as 1500K can be fitted with a low scaling factor, r\ ~0.05, although 
the resulting fit is not as good. 



Parameter 


value 


error 


units 


Component 1 


FWHM major 


20 


± 3 


mas 


FWHM minor 


6.4 


±1 


mas 


<P 


-42 


± 2 


degrees 


T 


800* 


± 150 


K 


r 


1.9 


± 0.5 




V 


0.25* 


± 0.07 




Component 2 


FWHM major 


56.5 


± 5 


mas 


FWHM minor 


42.4 


+5 


mas 


<P 





+70 

-20 


degrees 


T 


290 


± 10 


K 


T 


0.42 


± 0.2 




'1 


0.64 


± 0.15 
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Figure 7. Maximum entropy reconstruction at 8/im. Image size is 30 X 30 
pixels, with lpix=2 mas. Colour scale is linear. The extended blobs are arti- 
facts caused by the low uv coverage. Gaussian fitting to this image measures 
it to be 7.7 / 21 mas in FWHM with PA=-46°, very close values to the grey 
body model results, and in agreement with the results of the one Gaussian 
fitting at 8/im. The 'dirty map' is shown in the bottom left corner. 



4.1.1 Symmetry issues 

As stated before, only phaseless and therefore symmetric models 
were applied to the data. The lack of an absolute phase also implies 
a loss of astrometric information. We cannot, for example, deter- 
mine the position of the mid-infrared emission in relation to other 
known components on similar scales for NGC 1068. In principle, 
we could constrain the relative position of our two model compo- 
nents by means of the I, m parameters in equation[2] If our two com- 
ponents are not symmetric with respect to each other, the resulting 
phase difference will also affect the sum of the Fourier amplitudes 
of both components, i.e. the observed correlated fluxes. In practice, 
we find that we do not have enough information on component 2 in 
order to determine its relative position. Models with offsets of up to 
a few tens of mas are virtually indistinguishable from models with 
zero offsets. In order to better constrain the relative positions, more 
observations on shorter baselines are needed, as component 2 is 
mostly over-resolved in our baseline sample. We do know from the 
differential phase that the mid-infrared emission is asymmetrical. 
However, we cannot tell whether this asymmetry is due to the ge- 
ometry of the dust clouds or the result of asymmetrical absorption 
in the torus. 



4.1.2 The infrared SED 

The model used here is specifically designed to explain the MIDI 
data in the wavelength range of 8-13yum. It is not a full radiative 
transfer model, and therefore cannot be used to predict or account 
for the entire infrared SED of NGC 1068. 



4.2 Maximum entropy imaging 

The fits of the two blackbodies presented in 34. H are not perfect fits 
to the data. Attempting to improve them, we explored different per- 
turbations in the Gaussian distribution and although some of them 



can make the fit better, they are still ad-hoc. Finally we decided 
to reconstruct an image using maximum entropy (ME) methods, 
which guarantees that the resulting image will be the most statisti- 
cally probable reconstruction given the information in the data. To 
our knowledge, this is the first time such a method has been used 
with infrared interferometric data. The ME reconstruction is based 
on the general algorithm by Skilling & Bryan (1984) and modi- 
fied to the specific needs of dealing with MIDI data. The code was 
tested on simulated data with similar S/N as the original data and 
using the same uv coverage of NGC 1068 (Fig. [TJ. In general, the 
code is able to accurately reproduce Gaussian brightness distribu- 
tions or similarly simple shapes without sharp edges. Attempting to 
reconstruct multiple component images (such as a binary) have not 
been successful for this uv coverage. We have chosen to reconstruct 
the image at 8/jm, since in that wavelength a single component fits 
the data well and the silicate absorption does not come into play. 
The resulting image is displayed and discussed in Fig. [7] The ME 
method is useful since it is not model-dependent and no prior as- 
sumptions about the source are needed to make use of it. Although 
the ME image details the basic properties of the source, and can 
provide a good starting point for any model, it does not provide a 
more detailed picture than parametric models which enable us to 
probe the wavelength dependence of the data as well as its bright- 
ness distribution. Nevertheless, the ME image does fit the data set 
at 8yum perfectly. Since the code, by definition, always prefers to 
get rid of complicated structures, we can gain insights into the de- 
viations from the pure Gaussian shape in our models by looking at 
how the code chooses to reduce the x 2 of the image. The central 
shape in the image fits well with component 1 from 34.11 while the 
extended blobs seen are needed to perfect the x 1 - We believe that 
the blobs correspond to extended emissions present in the source, 
while keeping in mind that the specific shape and positions of the 
blobs mostly reflect the low uv coverage of the data set. The contri- 
bution of the extended emission, as measured from the intensities 
in the image, is about 15% of the flux of component 1. 
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Figure 8. Results of one Gaussian fitting, (a) FWHM sizes for the major 
(red) and minor (black) axes, overlaid with r oc ifx (dashed); (b) Gaussian 
peak; (c) position angle and (d) normalised % 2 , for each wavelength. We 
stress that in this model each wavelength was fitted independently, and the 
X 2 in panel (d) is the minimized^ 2 found for each wavelength. The 'bumps' 
near 9.7/im are caused by the O3 atmospheric feature. 



fore larger radii. For optically thin dust in radiation equilibrium 
where T oc r~ 05 , and for blackbody emission where A scales as 1/T 
(Wien's law), we get r oc s/a as a crude estimate, which fits the 
trend seen in the fitting of the Gaussian's axis. The discontinuity at 
A = 9.5yum is due to the atmospheric 9.7yum Ozone feature and the 
intrinsic Silicate absorption. 

(ii) The Gaussian flux density (Fig. [Hp) is simply a response to 
the varying flux with wavelength. 

(iii) The Gaussian's position angle (Fig. [8};) is perhaps the most 
interesting of results. First, we fit the same position angle from 8 
to 9/um, followed by a jump of ~10 degrees, which gradually ro- 
tates back but does not return to its original angle. The location of 
the jump in the position angle is indicative of the wavelength at 
which the Silicate absorption feature takes effect. The gradual, lin- 
ear change in the position can be most easily interpreted as asym- 
metrical absorption, i.e. the absorption at 9 micron is more pro- 
nounced to the south-west of the Gaussian than at 12yum. In con- 
trast, a constant PA will indicate a uniform absorption screen is 
effectively present. 

(iv) Further, we can deduce the presence of another dust compo- 
nent from the difference in PA between 8 and 13/um. The^ 2 (Fig. 
[Hp") also supports the presence of an additional component, grad- 
ually becoming worse from lOpm (with increasing wavelengths). 
Although we cannot determine the geometry of this second com- 
ponent from this method, we can estimate that it consists of warm 
(~300K) dust, whose blackbody emission begins contributing to 
the flux at A >lOfim. 

In general, the single Gaussian per wavelength model fits the 
data well, and helps to reduce the complex behaviour of the corre- 
lated fluxes to a set of four parameters that vary with A. It is, how- 
ever, a geometrical model and does not relate the different wave- 
lengths in a physical way. We find it encouraging that the different 
parameters are continuous and behave in an 'orderly' fashion, al- 
though fitted independently. 



4.4 The properties of each component 

We now discuss in detail the properties of each component in the 
model, combining results from the different modelling methods. 



4.3 Channel by channel Gaussian fits 

MIDI uses a grism to disperse the two beams into 261 channels be- 
fore combining them. Per channel per time, we have a set of 16 cor- 
related fluxes, or Fourier amplitudes, one for each baseline. From 
these amplitudes a Gaussian brightness distribution corresponding 
to infrared emission from this channel's wavelength is fitted. The 
fitting is done directly in the uv plane, and takes into account only 
the correlated fluxes. In each wavelength we fitted a single Gaus- 
sian, motivated by the successful Gaussian fitting of J04 to data 
composed of two baselines. The fit results are summarized in Fig. 
[8j which plots the values of the model's parameters (i.e. the ma- 
jor and minor FWHM sizes, height and position angle of the fitted 
Gaussian) and the resulting normalised x 2 f° r each fit, determined 
by dividing the^ 2 by the number of degrees of freedom (16) minus 
the number of model parameters (4). The main results from this 
approach are: 

(i) The FWHM major and minor axis (Fig. [8^) lengths in- 
crease with wavelength. This is expected of centrally heated dust 
since longer wavelengths indicate cooler temperatures, and there- 



4.4.1 Component 1 

Component 1 is the dominant source of Mid-infrared emission on 
MIDI'S scales. It is resolved in all of our observations and its geo- 
metrical properties are well constrained to have a major axis of 20 
and a minor axis of 6.4 mas in FWHM with a PA=-42°. As Table 
|2]shows, and in contrast to the second component, the temperature 
and optical depth are not well constrained. Nevertheless we can es- 
tablish that component 1 is composed of hot dust (T=:800K), with 
a low scaling factor (n ~0.2). Although we cannot unambiguously 
determine the chemical composition of the dust, we observe that 
it differs from the profile of standard interstellar dust. Namely, the 
dust should begin absorbing towards 9yum and the optical depth 
should rise more sharply towards 10/jm, followed by a shallower 
decline up to 14yum. Figure [8}; indicates asymmetric absorption 
toward lower position angles, which may explain the deviations 
from a 'perfect' fit in both our models for the wavelength range 
A >l0/im. Clumpiness: Evidence for dumpiness comes from the 
small scaling factor, which is consistent with clumpy emission on 
unresolved scales that effectively reduces the surface brightness of 
the source independently of wavelength. Attempts to fit the data 
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with models where the scaling factor is close to unity have not 
been successful. Extinction in the mid-infrared towards our line of 
sight probably contributes to the low value of 77, but cannot fully 
account for it without invoking dumpiness. More direct evidence 
for dumpiness comes from the profile of the correlated flux of uv 
#4 at 8 pm, which is significantly less than our Gaussian model pre- 
dicts. This is the observation with the longest baseline and as such 
it is most sensitive to small scale structure such as clumps. Given 
a resolution of 5.8 mas at 8/im for that baseline, we set an upper 
limit of 0.4pc (5.8 mas) for the size of a clump. Alternatively, one 
can interpret the missing flux at that uv point as a smooth (rather 
than clumpy), non-Gaussian flux distribution at sub-parsec scales. 

4.4.2 Component 2 

Since most of our observations were taken with long baselines, 
component 2 is mostly over-resolved. Figure[4](filled curves) plots 
the contribution of this component to the correlated fluxes, illustrat- 
ing its presence in only five of the sixteen baselines observed, cor- 
responding to those shorter than sixty meters in telescope separa- 
tion. As a result, its geometrical properties are not well constrained, 
while the temperature prediction is confirmed by that found in 34.31 
As for the size of the flux distribution, component 2 is an extended 
Gaussian structure, with a minor axis of 42 mas and major axis of 
56 mas in FWHM. It's orientation is best fitted as north-south, with 
a large uncertainty. Component 2 marks the colder, extended part 
of the torus-like structure, and is (together with comp 1) the un- 
resolved source reported by IPoncelet et alj J2007h and others (see 
introduction). Its north-south elongation is common with elonga- 
tion of the mid-i nfrared tongue (Figure [5} of NGC 1068, which 
iBock etalj 12000) attribute to re-emission by dust and UV radia- 
tion concentrated in the ionization cone. 

4.4.3 The physical distinction between the components 

In our model the two components are separate and distinct, each 
with a fixed temperature. In reality, we expect the dust, on average, 
to have a smooth temperature and density distribution contained in 
one structure. In this sense, the two components are an abstraction. 
However, our findings here, and in particular the quality of the fit to 
the data, indicate that the two-component approximation is accurate 
regarding the brightness distribution of the dust, and this requires 
a steeper temperature gradient than the simple T oc A expected 
for centrally heated, optically thin dust. 

4.5 Summary of modelling results 

To summarise our main findings, the mid-infrared emission from 
the core of NGC 1068 within a beam size of 28 mas (~2 pc ) is 
dominated by a warm Gaussian-shaped structure of dust (comp 1) 
with a chemical composition unlike that of dust towards our own 
galactic centre. The geometrical properties of this component are 
independently well constrained in each of our methods of investi- 
gation. Component 1 is co-linear with, and has similar size to the 
H2O megamaser disk. It is tilted by ~ 45° to the radio jet, and per- 
haps also by a lesser amount (~ 15°) to the ionization cone. The dust 
temperature is ~800K, assuming a grey body model. Evidence for 
dumpiness is mostly indirect, and comes from the low grey body 
scaling factor and the deviation from a Gaussian fit for data ob- 
tained with the longest baseline. This single Gaussian component 
cannot account for the total flux as observed with an 8m telescope, 



as shown in Figure [2] It also does not fit the data as well for wave- 
lengths longer that 10/jm on our shortest baselines, suggesting the 
need for a second component (comp 2), which is extended so that 
it will become over-resolved with baselines longer than 60m. The 
filled curves in Figure[4] which plots the contribution of component 
2 to the correlated fluxes, demonstrates how small that contribution 
is. As a result, the geometrical properties of this component, in- 
cluding its Gaussian nature, are not well constrained. 

4.6 Comparison with previous MIDI studies of NGC 1068 

Two previous studies of the core of NGC 1068 have been made 
using MIDI. Both studies were based on the same set of data: two 
tn'-points with baseline lengths of 42 and 72 meters, obtained dur- 
ing the science demonstration time of the instrument. They were 
taken with the PRISM mode, which offers a low spectral resolution 
of A/SA ^30. The first study, by Jaffe et al. (J04), used a simpler 
version of the two Gaussian grey-body components used in this 
work (j j4.lt . and their findings are summarised in the introduction. 
In general, our findings here agree with those of J04, deviating from 
their results only in one parameter: the orientation of the extended 
component (component 2), which was then reported to be elon- 
gated east-west. This discrepancy is not surprising in light of our 
finding here that the extended component is over-resolved in ob- 
servations with baseline lengths longer than 60m. This leaves just 
one observation (in J04) where component 2 contributes signifi- 
cantly to the correlated flux, not enough to determine its orienta- 
tion. The FWHM sizes of the extended component in Jaffe et al. 
(30 and 49 mas) are smaller than the sizes reported here, and so the 
component was able to contribute to the correlated flux obtained 
with the longer baseline (72m), and an estimation o f its orientation 
and ax is ratio could be made. The second study, bv lPoncelet et al.l 
(2006) (hereafter P06), used a different approach to modelling the 
same data set. Instead of Gaussian components, P06 used a com- 
bination of two uniform, circular disk components and visibilities 
rather than correlated fluxes as the quantity to be fitted. The derived 
angular sizes and temperatures of the two disk components are ~35 
and 83 mas, and ~361K and 226K, respectively. Both approaches 
provided good fits to the data from the two short baselines. This is 
because the differences between a Gaussian model and a uniform 
disk is small for short baselines. With the addition of data taken at 
longer baselines, we can now determine that the circular disk model 
of P06 is inconsistent with the MIDI data presented here. This is 
mostly due to the sharp edge of the circular disk model: the visi- 
bility (and the correlated flux) of a circular disk is proportional to 
the absolute value of the J-Bessel function of the first order, which 
goes down to zero at specific baselines depending on the angular 
size of the disk and the baselines length in units of the observed 
wavelength. For the disk sizes of the P06 model, these zeroes fall 
between baseline lengths of 60m (at 8/mi) to 100m (at 12.8/jm), a 
range which is well sampled by our observations. Yet, no trace of 
the drop in the correlated flux predicted by the uniform disk model 
is seen in any of our observations. The zeroes in the visibility func- 
tion appear for any flux distribution that possesses a sharp edge. 
Thus, our data is inconsistent with any such flux distribution. 



5 DISCUSSION 

NGC 1068 is considered as the prototypical Seyfert II galaxy. It 
is also the first Seyfert II with broad emission l ines seen in po- 
larized, scattered light dAntonucci & Milled !! 985 ) . These observa- 
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Figure 9. Picture summarising the multi-wavelength structures on parsec- 
scales in the nucleus of NGC 1068. The FWHM of the compact dust (com- 
ponent 1) is sketched in re d, centred around the H2O maser spots and 5GHz 
radio emission, both from Gallimore et al. (2004). Contour levels below the 
FWHM level have been removed to allow a better comparison between the 
ra dio and the mid-in frared. The ionization cones inferred from spectroscopy 
bv lDas et alJfeOOfj) are shown in yellow, and the HST [OIII ] (reduced in 
scale by a factor ~ 100) image contours (Evans et al. 1991) are shown in 
blue to indicate the visual orientation of the cone. 



tions proved that NGC 1068 harbours an obscured Seyfert I nu- 
cleus, and led to the AGN unification model, at least in its most 
simple form, uniting Seyfert I & II nuclei. This is the object for 
which the obscuring dusty torus was originally invented. The pre- 
vious sections have described the appearance of this structure at 
mid-infrared wavelengths and sub-parsec scales, and its geometric 
relation to structures seen at other wavelengths. 

The literature on the spectrum of NGC 1068, both observa- 
tional and theoretical, is very extensive. Many authors have mod- 
elled the production of the overall SED of the galaxy within the 
context of the unified theories. Without detailed structural obser- 
vations, however, they have found the problem underconstrained. 
The early interferometric work of J04 and W04 led modellers of 
the emitted spectrum in the direction of collections of optically 
thick inhom ogeneous ("clumpy") systems. Following the pioneer- 
ing wo rk of Nenkova et al. (2002) on clumpy models. Frlonig et all 
( 2006) and lSchartmann et al.l 020080 showed that three-dimensional 
radiative transfer models of clumpy tori were consistent with the 
low spatial resol ution SED data and the J04 interferometric data of 
NGC 1068. The iHonig eTal] j2006|) model of a torus with an in- 
clination of i = 55° is also able to fit the K-band mea surements of 
W04 t o a factor of 2 in correlated flux. In a later work (Honig et al. 
J2003)) the authors have revised the inclination angle of the torus 
to ( = 70°, this time providing excellent fits to both the SED and 
the previous interferometric N-band data, though still falling short 
of reproducing the K-band interferometric data. Mos t recently, a 
model with an edge-on torus (f = 90°) was presented (Honig et al. 
(2008)). This model was fitted to the SED and one of the interfero- 
metric data points of J04, reproducing a good fit to the former, and 
a fit correct to a factor of ~2 for the latter. Their most important 
conclusion is that the random cloud arrangement has a significant 
effect on the SED and images, providing remarkably different SEDs 
for different cloud arrangements and for the same set of large-scale 
torus parameters. This accounts for the difficulty in determining the 
inclination angle. 



The observational material presented here is considerably 
more detailed than that which was available to the authors cited 
above, and the results of their simulations of the production and 
redistribution of the radiative energy are generally not presented in 
a form that can be directly compared to our measured correlated 
fluxes as a function of wavelength and baseline and position an- 
gle. Hence, while these models seem globally consistent with our 
present data, we refrain from considering whether any of them do 
or do not fit it in detail. Rather we focus on the physical struc- 
tures we observe and their relations to the surrounding AGN com- 
ponents. We also consider the similarities and differences of the 
results f or NGC 1068 with t hose of the nearby Circinus Seyfert II 
galaxy, dTristram et alj|20"07r) . and compare the implications of our 
N-band data with those measured on the same galaxy by W04 in 
the K-band. 



5.1 The inclination and thickness of the dust 

In this section we argue that the nuclear dust structure in NGC 1068 
is circum-nuclear and inclined edge-on, and therefore its observed 
size and the 3:1 ratio of its major/minor axis is indicative of the 
thickness of its inner part rather than a projection effect of an 
inclined disk. Our main argument here is that the dust and the 
maser disk are co-spatial, and since the masers are seen edge-on, 
so is the dust. The possibility of detecting the hot dust through the 
cold dust that surrounds it has been demonstrated by several recent 
three dimensional radiat ive transfer calculation of tori: 
ISchartmann et al.l ( 120051) find that "in a homogeneous dust distri- 
bution the observed mid-infrared emission is dominated by the 
inner funnel of the torus, even when obser ving along the equ atorial 
plane". Similar results are also found by Honig et al. d2006h and 
discussed in more detail in section |4~6l We show that this scenario 
of an edge-on geometrically thick torus not only agrees with 3D 
radiative transfer models, but also is consistent with other observed 
structures on milliarcseconds scales: 

(a) The position angle of component 1 is virtually the same as the 
PA of the western maser spots. This can hardly be a coincidence. 
The theoretical relationship between water maser excitation and 
warm molecular dust has been well established (Neufeld et al. 
1994, for example), and the masers are indeed expected to be 
embedded in warm (>600K) molecular dust. 

(b) Apart from the position angle, the inner radius of the maser 
disk i s approximately the radius at which dust sublimes jGreenhilll 
1998). Dust reverberation modelling for the maser disk (Gallimore 
et al 2001, Fig. 9) has established the geometry of the masers, with 
those closest to the nucleus outlining a ring 0.6 pc in radius. The 
diameter of the ring (1.2 pc) is similar to the FWHM size of the 
major axis of our component 1 (1.35 pc), supporting the conclusion 
that the masers and the hot dust component are co-spatial. 



5.2 The orientation with respect to the jet 

Although there seems to be little or no correlation between the rel- 
ative angle of je ts and the dust res iding in the galactic disk for 
Seyfert galaxies fanney et al.ll2000T) . the naive expectation is that 
the dust directly orbiting the black hole and the accretion disk is 
perpendicular to the jet. It is the common belief that the orienta- 
tion of the jet is closely coupled to the spin axis of the black hole, 
which in turn is affected by the sum of all the angular momenta of 
the material accreted throughout its lifetime. ICapetti et alj i 19991) 
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estimate the lifetime of the jet as <1.5xl0 5 years, a much shorter 
time compared to the lifetime of the black hole. The current phase 
of the AGN is not likely to have affected the spin axis of the black 
hole, and it is not surprising that the dust and the maser disk are not 
aligned with the jet. Still, it is interesting to consider the origin of 
the orientation of the dust component and the maser disk, since it 
seems related neither to the jet nor to the galactic plane. We now 
consider two mechanisms: 

self warping. The sub-Keplerian velocities of the masers have been 
interpreted as evidence that the disk is superma ssive and will warp 
due to self-irradiation fcodato & Bertirj)2003h . This mechanism, 
however, only applies to thin structures and therefore cannot ac- 
count for the orientation of the dust. 

misaligned inflow. In this scenario the matter falls misaligned with 
respect to the current orientation of the jet. This can be the result 
of either matter infalling from outside the plane of the galaxy, or 
due to the influence of a torque, perhaps from a non-spherical nu- 
clear star cluster. The presence of such a potential will be most 
likely due to a past minor m erger, for which some evidence exists 
dGarcia-Lorenzo et alJI 19971) . Observations so far have determined 
that a dynamical mass of 6.5 x 10 8 M o is present inside a 50 pc 
radius dThatte et al Jl f 9971) . but due to obscuration by dust there is 
no information on the gravitational potential at s maller scales. This 
mass is about 12 times the mass of the black hole dLodato & Bertirj 
2003), and at a distance a few parsecs the black hole's potential still 
dominates over the potential of the star cluster. The influence of a 
star cluster may also explain the non-Keplarian velocity distribu- 
tion of the masers. Finally, the generally linear, stable shape of the 
jet excludes non-continuous capture of material into the black hole 
such as individual stars or individual dust clouds, which may take 
an arbitrary route. To conclude, the identical orientation of both 
the thin maser disk and the thick dusty torus strongly indicate that 
it is the angular momentum of the infalling matter combined with 
the gravitational potential at the nucleus that are responsible for the 
current orientation of the maser (i.e. accretion) disk and the inner 
part of the dusty torus that surrounds it. When discussing the ori- 
entation of the central dust component and the orientations of other 
components re l ated t o the nucleus, it is interesting to point out that 
Gallia no"et alj ([2003) have proposed to explain the asymmetrical 
appearance of the IT brightness distribution (on arcsecond scales) 
by suggesting that the Compton thick material located in the very 
central parts of the AGN is tilted with respect to the large scale 
molecular gas distribution with an orientation of 30 degrees, which 
is comparable to our findings. 

5.3 Relation between the torus and the ionization cone 

According to the AGN unification scheme, it is the obscuring torus 
that is responsible for the conical shape of the ionized gas in the 
narrow line region since it blocks the UV radiation from ionizing 
the material outside the cone. The cone's opening angle and orien- 
tation are then related to the orientation and geometry of the torus, 
and in particular to the arrangement of the clouds in its inner edge. 
Now that the inner part of the torus has been resolved, we can com- 
pare the two quantities and relate them to the predictions of the uni- 
fication model. The position angle of the dust is well constrained 
to be -42° (see Table |2j. We then expect the ionization cone to 
be centred perpendicularly (i.e. at PA=48°). The opening angle in 
this scenario depends mainly on the inner geometry of the torus, as 
each cloud in itself is optically thick to UV and X-rays, and trivial 
geometrical considerations limit the cone's angular extent to within 
±90° of the torus axis (i.e. -42<PA<138). 



The ionization cone, as seen in HST (Hubble Space Telescope) 
images, is centred around PA=10°, with a small opening angle of 
~ 45° (see Evans et al. 1991, for example). As we can see in Figure 
[9] its axis of symmetry is not aligned with the axis of component 
1, as it lies about 40 degrees further to the north. Yet, due to its 
small opening angle, the cone's angular extent is still within the 
bounds listed above. This discrepancy may be resolved by looking 
at the kinematics of the outflowing gas. lDas et alj J20061) have ap- 
plied bi-conical outflow models to high-resolution long-slit spectra 
of the narrow-line region obtained with the Space Telescope Imag- 
ing Spectrograph aboard the HST. These models require the cone to 
have a wider opening angle, ~ 80°, centred around PA=30°, reject- 
ing models with PA=10° (Das, private comm.), and in agreement 
with the orientation of our component 1. Figure [9] plots the HST 
image on top of the kinematic model. The northern edge of the 
cone coincides in both the model and the image, and it appears as if 
only half of the bi-conical outflow described by the model is seen in 
the HST image. The UV radiation which ionizes the gas in the cone 
is partially blocked by the obscuring dust closest to the AGN and 
is therefore sensitive to the configuration of the obscuring clouds. 
It is not unreasonable to suggest that some of the UV radiation is 
blocked by a cloud of dust so that only a part of the gas in the cone 
is ionize, as seen in the HST images. Assuming this, the orienta- 
tion of the dust component and the ionization cone roughly match 
our expectations and give further support to our interpretation of 
component 1 being the inner hot funnel of the torus-like obscur- 
ing object. This implies irregularities in the dust distribution and 
supports the idea of a clumpy dust distribution. 

If this scenario holds then we expect a difference between the 
appearance of the cone in the optical and in the infrared. The in- 
frared cone (i.e. narrow line region - NLR) should appear more 
symmetrical with respect to the dust since its appearance does 
not depend on the io nization of its material. The 12.8m image of 
iGalliano et all {2005) shows extended emission that is quite sym- 
metrical to the north and the south of the central engine, although 
the northern emission is considerably more pronounced. The emis- 
sion is centred on PA=~27 degrees, which is more symmetrical 
with respect to the dust, wit h a relatively small opening a ngle of 
-50 degrees. Alternatively, IPoncelet. Sol & Perrinl (|2008j) study 
the gas dynamic inside the NLR by looking at several mid-infrared 
emission lines by extracting spectra along two slit positions. They 
find the infrared cone to be centred on PA = 10° with an open- 
ing angle of ~80°. It seems then that in the infrared there is also 
disagreement between different methods of investigation. A case 
against this scenario of partial obscuration comes from the well de- 
fined shape of the cone over a scale of a few hundred pc, along with 
the linear shape of the radio jet both suggest that these irregularities 
in the dust distribution are stable over ~1500 years. 

To conclude, observations at infrared vs. optical wavelengths, 
and emission line analysis vs. imaging all disagree on the properties 
of the ionization cone / NLR. Further investigation is needed in 
order to clarify this puzzling aspect. 

5.4 The dust and the radio emission 

Radio observations at 5 and 8.4GHz of the nucleus of NGC 1068 
show a parsec-sized component (hereafter: 'radio component'), 
most likely indic ating free-free emissi on from the X-ray irradiated 
accretion disk iGalli more et alj|2004l) . If so, we expect the radio 
component to have similar geometrical characteristics to our com- 
ponent 1. Figure [9] shows the 5GHZ contours (up to the FWHM 
level) against the FWHM of component 1. The radio component is 
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'thicker' with a majonminor ratio of ~ 2 : 1, as opposed to 3 : 1 
for component 1 with its PA of -57 degrees. The two components 
are indeed similar in size and orientation. 



5.5 Comparison with Circinus 

It is interesting to compare our findings with those for the nucleus 
of the Circinus galaxy, the only other active galactic nucleus for 
which the central mid-infrared source was studied in similar de- 
tail. The Circinus galaxy is similarly regarded as one of the typical 
representatives of a Seyfert II galaxy population. It shows many 
of the characteristic attributes of Seyfert II sources, s uch as nar- 
row and b road emission lines dOliva et al. an ioniza- 
tion cone dVeilleux & Bland-Hawthornf 1997h, as well as a .4 pc 
edge-on, warped disk of HtO masers jGreenhilletai]|2003h . At 
about 4Mpc, the Circinus galaxy is located significantly closer to 
us than NGC 10 68. However, due to i ts lower nuclear luminosity, 
L acc = 10 10 L G dTristram et al.l 12007). it appears slightly fainter 
in the mid-infrared than NGC 1068. Using an approach to anal- 
yse the MIDI data s imilar to the one used by us for NGC 1068, 
iTristram et al.l d2007t) have found the dust in the nucleus of the 
Circinus galaxy to be distributed in two components: (1) a dense, 
warm (~330K) component of 0.4 pc size and (2) a slightly cooler 
(~300K), geometrically thick torus component with a size of 2.0 pc. 
The model fits to the data were improved by introducing clumpy- 
like perturbations in the flux distribution, mainly of the larger com- 
ponent, thereby providing arguably the first direct evidence of a 
clumpy structure in an AGN torus. 

As for NGC 1068, the compact component coincides with the 
nuclear water maser disk in orientation and size and there is evi- 
dence for a clumpy or filamentary distribution of the dust. Hence 
there seems to be some similarity between the dust distributions in 
these two Seyfert galaxies. Especially striking is the agreement be- 
tween the sizes when considering that the torus size is expected to 
scale with the square root of the luminosity of the central energy 
source, and the identical orientations of the dust and the masers. 
However, when the properties of the dust distributions are com- 
pared in more detail, there are several aspects in which the dust 
distribution in NGC 1068 deviates from that of the Circinus nu- 
cleus significantly. In the first place, while the temperatures of the 
extended components in the two galaxies are comparable, the com- 
pact component in the Circinus galaxy is significantly cooler (T = 
330 K) than the one in NGC 1068 (T = 800 K). The shallow tem- 
perature decrease of the dust in the Circinus gal axy was attributed 
to a high degree of dumpiness in the torus ( Sch artmann et al.l 
2008). In fact, this torus seems to be lacking a significant amount of 
truly hot dust. And indeed, we are observing different physical dust 
structures: in NGC 1068 the compact component is interpreted as 
the dust in the inner funnel of the torus, which is heated up to near 
sublimation temperature. In the Circinus galaxy, on the other hand, 
the compact component is interpreted as a disk-like dust structure 
in the centre (or 'midplane') of the otherwise geometrically thick 
torus. This structure is co-spatial with the rotating disk of maser 
emitters. As a second difference, in the Circinus nucleus, there is 
less silicate absorption in the compact component than in the ex- 
tended component. This was interpreted to be due to the silicate 
feature in emission, as is expected from optically thin dust in the 
inner regions of the torus. In NGC 1068 the inverse is the case, 
the absorption feature deepens for the compact component. This, 
and the unusual spectral characteristics of the dust in NGC 1068, 
as opposed to the dust in the Circinus galaxy that was found to fit 
that of standard galactic dust, suggest the reprocessing of dust in 



NGC 1068, in analogy to the evolution of dust in protoplanetary 
disks. Several underlying properties can have a major influence on 
the observed temperatures and absorption depths: among them are 
the total luminosity of the central energy source, the exact inclina- 
tion angle of the torus, and the volume filling factor of the torus. 
Proper radiative transfer calculations for the dust distributions in 
these two Seyfert galaxies are needed to bring more light to this 
issue. The relationship between the inner dust component and the 
maser disk found in both galaxies is the most striking of the fea- 
tures both objects have in common, and it may be that this is in- 
deed the generic case in all Seyfert nuclei that also possess maser 
disks. To conclude, the existence of a thick configuration of ob- 
scuring dust was confirmed in both Circinus and NGC 1068. Both 
dust configuration are found to fit a two component Gaussian model 
with differences mostly in the temperature of the components, and 
both components are similarly related to the water maser disk and 
the ionization cone. 



5.6 Comparison with other infrared interferometric 
measurements 

This galaxy has been observed previously in th e infrared with in - 
terferometric resolution in the N-band by us d-Taffe et all d2004l) ) 
and in the K-band by Wittkow ski et al. with the VINCI instru- 
ment, ( Wittkows ki et al J ( 120041) W04), and in speckle mode by 
Weig eltet al.l d2004r) . These last observations were restricted to ef- 
fective baselines shorter than 6 meters by the telescope used, cor- 
responding to ~ 24 meters at our shortest baseline at 8 /um wave- 
length. Thus they are primarily sensitive to structures larger than 
any observed here, and we will not discuss them further. Here we 
summarize the findings of W04 and discuss possible nuclear struc- 
tures that are consistent with both their observations and those re- 
ported here. 

W04 found a squared visibility of 0. 16±0.04 at a projected 
baseline length of 45.8 meters and position angle of 44.5 degrees, 
corresponding to a resolution of ~5 mas. The authors favour a 
multi-component model where 52 mjy originates in an unresolved 
source of size <5 mas and 75 mjy arises on scales of the order of 
40 mas or larger. Again this larger component will not be discussed 
here. The smaller emission was attributed to thermal emission from 
hot (1000-1500 K) dust at the inner cavity of the torus, with a pos- 
sible contribution of direct plus scattered light from the central en- 
gine. The resolution achieved by VINCI at 2.2 jim is only ~ 30% 
smaller than that we reach on our longest baseline, 129 meters at 
our shortest wavelength, 8 /urn, thus in spatial terms the measure- 
ments can be compared directly. 

The interpretation of these measurements is clouded by sev- 
eral uncertainties. First, the interferometric observations provide no 
direct positional information. Thus we do not know the position of 
the small K-band component relative to the galactic nucleus or rel- 
ative to the extended features seen with MIDI, except that it must 
lie within the 56 mas field of view of the VINCI fibers, presum- 
ably centered on the nucleus. Second, the extinction in front of the 
K-band component is uncertain. If this component is seen through 
the same extinction as our Component 1 (Table 2), then the K-band 
absorption may be estimated by scaling the de pth of the silicate 
absorp tion by a standard extinction curve (c.f. ISchartmann et al.l 
(2005), Fig. 3), yielding an optical depth of t k ~ 3 but with large 
uncertainties, particularly if the dust is of non-standard composi- 
tion. If the K-band component is not at the nucleus, e. g. represents 
a clum p of warm dust clouds as suggested by W04 and lHonig et al.l 
(2008), then perhaps our line of sight to it does not pass through 
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the main body of dust absorption, and the extinction could be sig- 
nificantly smaller. Third, as we have seen in Section 4, the short 
wavelength, long baseline MIDI results can be well modelled by a 
smooth gaussian of size 20x 6 mas, with no evidence of unresolved 
emission corresponding to that found by VINCI. However, the low 
scaling factor of Component 1, a = 0.25 (Table 2), points to the 
possible existence of substructure on size scales below the resolu- 
tion of MIDI. For the 8 pm flux of any of the sub-components we 
can then only set an upper limit of ~ 750 mjy corresponding to the 
correlated flux seen on a 129 meter baseline (Figure 3, panel 1). 
Assuming again a standard extinction law, its intrinsic 8 pm flux 
would be < 1.3 Jy. 

Under these circumstances, we consider only two possibili- 
ties for the K-band source where reasonable additional assumptions 
can be made concerning the above uncertainties: either emission 
from the accretion disk around the central black hole, or from a 
small body of hot dust, as discussed by W04 and also proposed by 
Iffimig et al.H2008T) . 

5.6.1 Central Accretion Disk 

In this case we assume that the source lies behind the absorption 
seen in N-band and assume that for the intrinsic flux from the ac- 
cretion disk S,,(2.2yu) = 0.05 x e 3 = 1 Jy and 5,,(8yu) < 1.3 Jy. This 
would correspond to a v~°- 2 spectrum, not unreasonable for an ac- 
cretion disk. In addition, this accretion disk flux is consistent with 
the total accretion disk luminosity: From the total infrared emission 
j vS v d In v ^ 2.5 10 15 Jy Hz (estimated from the data assembled in 
the NASA/IPAC ExtraGalactic Database) one infers the bolometric 
emission of the disk to be vS v (peak) ~ 7 x 10 15 Jy Hz (assuming 
that the UV luminosity exceeds the IR luminosity by a factor ~ 3, 
because light emitted out of the torus plane is not reprocessed into 
the infrared). At a distance of 14 Mpc this corresponds to ~ 1.5 10 45 
erg s -1 , sim ilar to other estimates of the bolometric luminosity of 
NGC 1068. ISchartmann et all d2005l) discuss likely SEDs for cen- 
tral accretion disks, based on observations of Type I AGNs (e.g. 
Walter et al., 1993, Zheng et al., 1997) and more modern theoretical 
models where observations are lacking. In most cases the spectrum 
is relatively flat in the optical-IR region, -0.5 < a v < +0.3 (where 
S y oc v Bl ) and breaks sharply down shortwards of the Lyman edge 
(A-p ea k ~ O.lyum i.e. v peak ~ 3 10 15 Hz) With this value of v peak and 
the bolometric luminosity we find S v (peak) ~ 2 Jy. The relatively 
flat spectrum thus inferred between v peak and 2.2 pm is consistent 
with the acceptable range of a v . 

We conclude, therefore, that the fluxes and limits measured 
with the VLTI at 2 and 8 pm are consistent with an accretion disk 
spectrum. In fact, if the emission spectrum of the accretion disk 
and the (relatively low) foreground extinction value are correct we 
expect to see the accretion disk in the K-band measurements. Q 

5.6.2 Dust Cloud Emission 

Several authors (Nenkova et al., 2002, Honig 2006, Dullemond et 
al 2005) have proposed that AGN dust tori are clumped and Honig 
(2008) suggest that the small component of the VINCI observa- 
tions is a warm dust cloud or compact cluster of such clouds. To 

1 This is at variance with Honig et al. 12008), who on the basis of a sim- 
ple theoretical v 1 ' 3 accretion disk spectrum peaking at v ~ 10 17 Hz, for 
which little observational evidence exists, argued for an undetectably low 
accretion disk flux at A > Ifim. 



evaluate this suggestion in the light of our observations requires an 
additional assumption about the K-band extinction. If, as above, we 
assume the clouds to be near the nucleus this would to be of order 
> 3 mag, the intrinsic K-band flux is > 1 Jy and the 8pm — > 2.2pm 
spectral index is a v < -0.2. For black-body emission this implies a 
color temperature T co \ or > 1300 K, which is rather warm for dust, 
but not excluded, particularly if the dust is more refractory than 
standard silicates. If, in fact, the dust clouds are in the immediate 
vicinity of the nucleus, i.e. r < 2.5 mas, their equilibrium tempera- 
ture derived from the above luminosity would be much higher than 
this, making dust emission very unlikely. However, as discussed at 
the beginning of this section, there is no direct evidence that the 
emission feature is physically so close to the nucleus. If the unre- 
solved K-band source is not centered on the nucleus but above the 
main body of dust clouds, or seen through a hole in these clouds, 
then its intrinsic K-band flux could be as low as the observed 50 
mjy, in which case the N-band upper limit provides no useful in- 
formation on its properties. 

Thus we conclude that the K- and N- band interferometric data 
can be straightforwardly be explained as emission from a hot cen- 
tral accretion disk if reasonable assumptions are made on the fore- 
ground emission and the accretion disk spectrum. The data can also 
be explained as emission from a group of hot dust clouds, but then it 
is unlikely that these clouds are centered around the central source. 
Additionally in this case, the extinction towards the core must be 
sufficient (i.e. r K » 3) to suppress the expected K-band emission 
from the nuclear accretion disk below the observed 50 mJy level. 



6 CONCLUSIONS 

In this paper we present interferometric observations of the nucleus 
of NGC 1068, using MIDI at the VLT. Extensive uv coverage of 
sixteen baselines with a maximal resolution of 7 mas has allowed 
us to analyse the mid-infrared (8-13/jm) emission from the obscur- 
ing torus in great detail. We find the measurements consistent with 
emission from a compact (0.45 x 1.35 pc) component with a Gaus- 
sian flux distribution that is composed of hot (~800K) silicate dust 
with an unusual absorption profile, which we identify as a funnel 
of hot dust associated with the obscuring torus. The emission is 
co-linear, and likely co-spatial, with the well studied H 2 mega- 
maser disk, and therefore tilted by 45° with respect to the radio jet. 
A second, more extended (3x4 pc) component of warm (~300K) 
silicate dust is weakly detected, which we identify as 'body' of the 
torus. This second component is mostly over-resolved and its prop- 
erties are not well constrained. We discuss the physical origin of the 
emission with respect to the torus, the masers, the ionization cone, 
and the radio jet. A direct image of the source at 8 pm, obtained 
with maximum entropy image reconstruction, is presented as well. 
Since the obscuring torus is a crucial component in the accepted 
picture of an AGN, resolving the structure of the torus illuminates 
many aspects of the AGN picture. We show that in many aspects 
the nucleus of NGC 1068 is irregular: the orientation of the dust is 
tilted with respect to the jet; the PA of the visible ionization cone 
does not match the PA of the inner dust funnel; and the chemical 
composition of the dust in the torus is different than dust observed 
in the Milky Way. 
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